Enhancing synchronization in complex networks of coupled phase oscillators 
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By a model of coupled phase oscillators, we show analytically how synchronization in non-identical complex 
networks can be enhanced by introducing a proper gradient into the couplings. It is found that, by pointing 
the gradient from the large-degree to the small-degree nodes on each link, increase of the gradient strength will 
bring forward the onset of network synchronization monotonically, and, with the same gradient strength, hetero- 
geneous networks are more synchronizable than homogeneous networks. The findings are tested by extensive 
simulations and good agreement are found. 
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Synchronization in complex networks has been a topic of 
arising interest in recent years, mainly due to its implications 
to the practical processes observed in biological and neural 
systems uJ, |2]- While most of the studies are focusing on 
the phenomenon of complete synchronization in networks of 
identical node dynamics [3], there are also interests in the col- 
lective behaviors in non-identical networks |4, 5]. The model 
of non-identical network is more representative to the realistic 
situations and, to analyze its dynamical properties, requires 
some special mathematical methods |4, 5]. Different to the 
studies in identical networks, in non-identical networks peo- 
ple are usually interested in the onset of the system coherence, 
i.e. the critical coupling from where the systems transits from 
the incoherent to coherent states |5, 6]. For general node dy- 
namics, the critical coupling can be estimated by the method 
presented in Ref. |4], based on the information of the node dy- 
namics and the largest eigenvalue of the adjacency matrix of 
the network. However, for the special case of coupled phase 
oscillators, the onset of network synchronization could be de- 
scribed more accurately by some other approaches. For in- 
stance, it has been shown that the critical coupling character- 
izing the onset can be efficiently predicted based on only the 
information of network degree distribution, i.e. the mean-field 
(MF) approach in Ref. [5]. In all these studies, the network 
couplings are considered as of uniform strength, i.e. the un- 
weighted networks. 

Noticing that couplings in realistic networks are usually di- 
rected and weighted and, in many cases, the direction and 
weight of the couplings are determined by a scalar field [7], it 
is thus natural to extend the study of non-identical network to 
the weighted case. For identical networks, it is shown that the 
synchronizability of a complex network can be significantly 
improved by introducing gradient into the couplings 1 8ll9l [l0tl . 
So far these findings are obtained from identical networks and 
refeiTing to the transition of global network synchronization. 
In this paper, we are going to study the effects of coupling 
gradient on the onset of system coherence in non-identical net- 
works, and explore their dependence to the network topology. 
The former study will extend the MF approach of Ref. |5] 



to the situation of weighted network and the later stay could 
provide theoretical support to the findings of synchronization 
pathes in Ref. |11]. 

We consider network of N coupled phase oscillators of the 
following form (the generalized Kuramoto model ISJ) 
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with 0n and ojn the phase and natural frequency of oscillator 
n respectively, e is the overall coupling strength, and Cnm is 
an element of the coupling matrix C. In general, the matrix C 
is asymmetrical and the frequency aj„ follows some probabil- 
ity distribution p{uj). For the purpose of theoretical tractabil- 
ity, we assume that the network is densely connected and has 
a large size. Defining the global order parameter as r = 
Eti-n/EtiC, withr„e''/'" ^ El^iCn^ie''-), 
the local order parameter and d™ = J2m=i ^nm the total 
incoming couplings of n, then the onset of the network syn- 
chronization is characterized by the critical coupling strength 
Ec at which r starts to increase from 0. By the approaches of 
Ref. pO, we are able to obtain a similar equation for r in the 
region of e > £c 
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is the total outgoing couplings de- 

a\ — 2/[7rg(0)] and 0:2 = —t^g (0)ai/16 



with d°f = 
parting from n 

are two parameters determined by the first-order and second- 
order approximations of the frequency distribution piuS), re- 
spectively. The critical coupling is given by the following 
equation 
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with (. . .} denotes the system average. Please note that in our 
weighted model, the total incoming couplings d'" and the total 



outgoing couplings (i°"* of each node are real and, in general, 
unequal. Our main task is to investigate how the distributions 
of d™ and d°"* will affect the onset of network synchroniza- 
tion. 

We start by considering an unweighted, symmetrical net- 
work described by adjacency A = anm, with anm = 1 
if nodes n and m are connected, a„m = otherwise, and 



with 



0. The degree of node n is A:„ = ^ 



N 

m=l "n™ 



. To 

introduce gradient into the couplings, we transform matrix A 
as follows. For each pair of connected nodes n and m in the 
network, we deduce an amount g from a„m (the coupling that 
n receives from m) and add it to a„„i (the coupling that m 
receives from n). In doing this, the total couplings between n 
and m is keeping unchanged. Therefore a coupling gradient is 
generated which is pointing from node n to node m. Denot- 
ing the resulted matrix as 5*. The coupling matrix C is then 

defined as Cmn = knSnm/J2i=i ^nj for the non-diagonal 
elements, and C„„ — /c„ for the diagonal elements. The cou- 
pling gradient from n to m thus is AC„m = Cmn — Cnm ~ 

kmSmn/ A^j=i Smj~'knSn7n/ z^j=i ^nj- In realistic systcms, 
the direction and weight of each gradient are generally de- 
termined by a unified scalar field which, in the sense of net- 
work s^chronization, is usually defined on the node degree 
M> S [lOtl . Without losing generality, we make the gradient 
point from larger-degree to smaller-degree nodes on each link 
(the inverse case can be achieved hy g < 0). 

Now we discuss how the change of the gradient parameter 
g will affect the network synchronization. Noticing that in Eq. 
(O the value of a is independent of g and, by the definition of 
C, we always have (d™) — (k), which is also independent 
of g. Therefore the introduction of gradient will only affect 
the value of (^d™d°'^*y Rearranging the node index by a de- 
scending order of their degrees, i.e. ki > fc2 . . . > k^, then 
the outgoing couplings of n can be divided into two groups. 
Neighbors of node index m < n have element s„„i = 1 — g 
in matrix S and Cmn < 1 in matrix C (gradient points to n), 
while for nodes of index we have Smn = 1 + 5 in matrix S 
and Cmn > 1 in matrix C (gradient points to m). By this 
partition, the total outgoing couplings d""* is approximated as 
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with Pi<„ (Piyn) the probability for a randomly chosen 
node to have degree larger (smaller) than node n. Vti = 
■^ 'Yl,n=i Sij is the normalizing factor defined on node. In cal- 
culating Hi, again, we can divide the neighbors of i into two 
groups. Nodes of index j < i have Sy = 1 + g and nodes of 
index j > i have Sij = 1 — g. Based on this partition, we 
write 

a-l + .g(P,<, -P,>,), (5) 

with P the same definition as that of Eq. ^. For heteroge- 
neous networks of degree distribution P{k) — Ck~'^ and , 

we have r!, = 1 + ^ [2k]-'^ - k^-J - k'^,;^}, with kma. 
and kmin denote the largest and smallest node degrees of the 
network, respectively. Inserting this into Eq. (|4|i, we obtain 

C"* =k^[F + G„] (6) 
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Finally we have 
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Eq. (|9|l is our main result which tells how the network syn- 
chronization {Ec) changes with the coupling gradient ig) and 
the network topology (7). 

From Eq. (|6]l we know that, in comparison with the un- 
weighted networks, the introduction of the coupling gradient 
changes only the weight H = F + G of the outgoing cou- 
plings on each node, while in this process the total strength 
of the outgoing couplings is keeping unchanged. That is to 
say, gradient makes the distribution of H change from an 
even form {H = 1 in unweighted network) to an uneven form 
(H = H{g, k) in weighted network). Physically, the term F 
can be understood as a summation of the symmetrical part of 
the couplings on each node, i.e. F ^^ J2i=i iiiin(C„i, Ci„), 
which only depends on parameter g and will be decreased as 
g is increased. In contrast, the term G is a joint function of 
g and /c„. While G increases with g, its exact value, how- 
ever, are strongly modified by the node degree: larger degree 
assumes larger G [Eq. ([8]l]. The joint effect of F and G will 
divide the nodes into two groups. Nodes of degrees larger than 
some critical value kc have weight H > 1, while nodes of de- 
grees smaller than kc have weight H > 1. The critical degree 
kc can be calculated from the equation of i? = 1. Under the 
assumption of kmax » kmin, we have 
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The uneven distribution of H can be further understood by 
considering its approximations at /c « kmax and k w 
kmi7i, which results in Hu^k,^^^ = ^ (1 + .g) In j^ and 
Hk^kr^^r, = ^ (1 " 5) In i±|. Clearly, we have Hk^k^,^ > 
HkKk^in- Since the sum of H over the network is fixed, 
i.e. X]i=i ^i — ^' '^he gradient effect thus can be roughly 
regarded as a shifting of weight H from smaller-degree to 
higher-degree nodes. 

For scale-free networks generated by the standard BA 
growth model (ill, we have kmax ~ kminN^i^^ . Inserting 
this relation into Eq. dSj we obtain 



H = F -In 
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which basically tells the following: for fixed gradient parame- 
ter g, increasing the homogeneity of the network, i.e. increas- 
ing exponent 7, will make the distribution of H more homo- 
geneous and, as a result, the network synchronization will be 
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FIG. 1: (Color online) For scale-free network of 1500 nodes, aver- 
age degree 400, and degree exponent 7 = 3, the variation of the 
squared order parameter r as a function of the coupling strength 
£ in the region of e G [0.5£c, 1.5£c] by using gradient parameters 
g — 0.5 (the left symbol curve), g — (the middle symbol curve), 
and g = —0.5 (the right symbol curve). Apparently, the onset point 
of synchronization is shifted to the small values as g increases. Each 
data is averaged over 10 network realizations. The three line curves 
are plotted according to Eq. ^, which predicts the behavior of r^ 
reasonably well in the region of e £ [ec,0.3£c]- In all the three 
cases, the numerical results of the critical couplings £c are in good 
agreements with the theoretical predictions calculated from Eq. ^. 



suppressed (i.e. the value of Ec will increase with 7). Eq. (fTTl i 
gives the dependence of network synchronization on network 
topology. 

Now the effect of coupling gradient and the effect of topol- 
ogy on the starting of synchronization in nonidentical net- 
works can be summarized as follows. The changes of the gra- 
dient strength g or the degree exponent 7 do not change the 
total coupling cost of the network, they will only redistribute 
the weights of the outgoing couplings at each node according 
to its degree information. By adding gradient, the outgoing 
couplings at the small-degree nodes (of degree k < kc) will 
be reduced by an amount and added to those of large-degree 
nodes (of degree k > kc). This will induce a heterogeneous 
distribution in H which in turn will decrease the threshold 
coupling £c (see Eq. Q). This enhancement of network syn- 
chronization, however, is modulated by the network topology. 
By increasing the degree exponent 7, the distribution of H 
tends to be homogeneous (i.e. H ^ 1) and, consequently, 
network synchronization is suppressed. These are the mecha- 
nisms governing the effects of gradient and topology on net- 
work synchronization. The above analysis shows that: 1) the 
synchronization of non-identical networks can be enhanced by 
coupling gradient; and 2) in comparison with homogeneous 
networks, heterogeneous networks take more advantages from 
the coupling gradient. 

We now provide the numerical results. The networks are 
generated by a generalized BA model [12], which is able to 
generate networks of varying degree exponent 7. The fre- 



FIG. 2: For a scale-free network of 5000 nodes, average degree 100, 
and degree exponent 7 = 8, the variation of the critical coupling 
strength Ec as a function of the gradient parameter g. The solid line 
represents the theoretical results predicted by Eq. Q. 



quency distribution is given by piu) = (3/4) (1 — w^) for 
— 1 < cj < 1 and p{uj) = otherwise. The initial phase 9 
of each oscillator is randomly chosen within range [0, 27r]. A 
transition time T = 100 is discarded, and the value of r^ is 
calculated over another period of T = 100. To show the gra- 
dient effects on network synchronization, we have calculated 
the variations of the squared order parameter r^ as a function 
of the coupling strength e for three different gradient parame- 
ters: g — 0, 0.5 and —0.5. (According to our definition, .g < 
means that gradient is pointing from smaller to larger nodes.) 
The results are plotted in Fig. [T] Clearly, the critical coupling 
strength Sc is shifted to small values as g is increased. The 
three lines plotted in Fig. [T] represents the theoretical results 
of Eq. (|2|i, which fit well with the numerical results in the 
neighboring region of the onset. More importantly, the posi- 
tion of the onset coupling Sc is predicted precisely by Eq. (O. 
(The precision of this predication is dependent on the size and 
connectivity of the network, larger and denser networks give 
better results.) 

To have a global picture on the gradient effect, we plot Fig. 
|2]the simulation result of the variation of Sc as a function of 
g. It is shown that, as g changes from —1 to 1, the value of ec 
is monotonically decreased. This process of synchronization 
enhancement is well captured by Eq. ([8]), especially in the 
region of g > 0. Since in our analysis we have assumed the 
network to be of very large size and of dense connectivity, 
the mismatch between the theoretical and numerical results in 
Fig. |2]is reasonable. 

Simulations have been also conducted on the dependence 
of Ec on 7. By the generalized BA model II12I1 . we vary the 
degree exponent 7 continuously from 3 to 25, while keeping 
the size and average degree of the network unchanged. As 
we have predicted [Eq. (fTTTil. in Fig. [3] it is found that, for 
each value of g, the critical coupling Ec will increase mono- 
tonically with the degree exponent 7. Specially, for the case 
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FIG. 3: (Color online) For a scale-free network of 5000 nodes, av- 
erage degree 100, the variation of the critical coupling strength Sc as 
a function of the degree exponent 7 under the gradient parameters 
g = 1x10^^ (the upper symbol curve) and g = 0.2 (the lower 
symbol curve). For both cases, Sc increases with 7. The solid lines 
are the theoretical results predicted by Eq. dllb . 



network. The amazing thing is that, for network of given de- 
gree distribution (not limited to the scale-free type), our the- 
ory tells how much improvement could the network benefits 
from a given gradient. From the findings, we are able to not 
only point out clearly the optimal configuration for synchro- 
nization, which happens when g = 1 [Fig. |2l, but also have a 
systematic understanding on the transition from unweighted 
to optimal network, and, more importantly, the underlying 
mechanisms that govern this transition. Secondly, although 
similar findings about the gradient effects had been discovered 
previously in the study of global synchronization of identical 
networks |18|,|9|,11P], our analysises, however, are focusing on 
the onset synchronization in non-identical networks. Another 
difference is, by adopting the generalized Kuramoto model, 
that we are able to show analytically how the coupling gradi- 
ent affects synchronization (see Eq. ^ and Fig. |2]) and what 
is the role of network topology in this process ((see Eq. ( fTTI ) 
and Fig. [3]). It is noticed that in Ref. [| 11]] the authors found 
numerically that the onset of synchronization in scale-free net- 
works happens in advance to that of homogeneous networks, 
which, according to the approximation of Eq. ( fTTI ). can be 
easily understood. 



of g = 1 X 10^^ in Fig. [3] the numerical results are in good 
agreements with the theoretical results of Eq. JTOl l. As g in- 
creases the mismatch between the theoretical and numerical 
results is enlarged, especially for networks of larger 7. Again, 
by increasing size and coupling density of the network, the 
mismatch can be alleviated. 

A few remarks are in order. Firstly, while our theory gives 
well approximations on the collective behavior in densely 
connected large networks, our findings about gradient effects 
of their dependence to network topologies are general for any 



In summary, we have studied the effects of coupUng gradi- 
ent on the onset of synchronization in nonidentical complex 
networks. It is found that: 1) network synchronization can be 
enhanced by introducing gradient into to the couplings; and 2) 
in terms of the onset of synchronization, heterogeneous net- 
works are more synchronizable than homogeneous networks. 
We hope these findings to be helpful in understanding the col- 
lective behaviors in realistic systems. 

YCL was also supported by AFOSR under Grants No. 
FA9550-06- 1-0024 and No. FA9550-07- 1-0045. 
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By a model of coupled phase oscillators, we show analytically how synchronization in non-identical complex 
networks can be enhanced by introducing a proper gradient into the couplings. It is found that, by pointing 
the gradient from the large-degree to the small-degree nodes on each link, increasing the gradient strength 
will bring forward the onset of network synchronization monotonically, and, under the same gradient strength, 
heterogeneous networks are more synchronizable than homogeneous networks. These findings are verified by 
extensive simulations. 

PACS numbers: 89.75.-k, 05.45.Xt 



Synchronization of complex networks has received many 
interests in recent years, mainly due to its important implica- 
tions to the processes in biological and neural systems |lll,|2]- 
While most of the studies are focusing on the phenomenon 
of complete synchronization in identical networks (the same 
dynamics for all nodes) yl], there are also interests on the col- 
lective behaviors of non-identical networks (node dynamics is 
non-identical) |4, JJ. Comparing to identical networks, non- 
identical networks are more representative of the realistic sit- 
uations and, to analyze their dynamical properties, a set of 
special mathematical treatments had been developed fl ISfl- 
Meanwhile, in non-identical networks people are more inter- 
ested with the onset of network synchronization, i.e. the criti- 
cal coupling where the system transits from incoherent to co- 
herent states |l5l|6|], instead of global synchronization studied 
in identical networks. In general, the onset coupling strength 
of a non-identical network can be estimated by the method 
proposed in Ref. |4], based on the information of node dy- 
namics and network topology. However, for coupled phase 
oscillators, this estimation can be greatly simpUfied and im- 
proved by some special approaches. For instance, by the 
mean-field (MF) approach of Ref. ^, the onset coupling 
strength can be predicted by only the information of degree 
distribution. In all these studies, the network couplings are 
considered as having the uniform strength, i.e. the networks 
are unweighted. 

Noticing that couplings in realistic networks are usually di- 
rected and weighted and, in many cases, the direction and 
weight of each coupling are determined by a scalar field of the 
network 1 7], it is natural to extend the studies of non-identical 
networks to the weighted cases. For identical network, it has 
been shown that the synchronizability of a network can be 
significantly improved by introducing some proper gradients 
into the couplings |8, 9, 10]. So far these findings are ob- 
tained from identical networks and referring to the transition 
of global network synchronization. In this paper, we are going 
to study the effects of coupling gradient on the onset of system 
coherence in non-identical networks, and explore their depen- 
dence to the network topology. The former study will extend 



the MF approach of Ref. |J5|l to the situation of weighted net- 
works, and the later study could provide insights to the prob- 
lem of synchronization pathes observed in Ref. lUll] . 

We consider network of N coupled phase oscillators of the 
following form (the generalized Kuramoto model ||5J) 
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with 9n and ojn the phase and natural frequency of oscillator 
n respectively, e is the overall coupling strength, and c,„„ is 
an element of the coupling matrix C. In general, the matrix C 
is asymmetrical and the frequency aj„ follows some probabil- 
ity distribution p{uj). For the purpose of theoretical tractabil- 
ity, we assume that the network is densely connected and has 
a large size. Defining the global order parameter as r = 
Ell rj Ell d^\ with r^e''^" ^ E™=i ^nm {e^'-), the 
local order parameter and d™ = Em=i Cnm the total in- 
coming coupling strength of node n, then the onset coupling 
strength Ec of network synchronization is defined as the point 
where r starts to increase from 0. By the similar approaches 
as used in Ref. ||5|], we are able to estimate the value of r in 
the region of e > Ec 
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In Eq. ^, d°"* = Em=i ^rnn denotes the total outgoing 
coupling strength of node n, ai = 2/[7Tp{0)] and a2 = 
—Tip (0)ai/16are two parameters determined by the first- 
order (p(0)) and second-order (p (0)) approximations of the 
frequency distribution, respectively. The critical coupling e^ 
reads 



ai 
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with (. . .) denotes a system average. Please note that in our 
weighted model, the total incoming coupling strength d™ and 



the total outgoing coupling strength d""* at each node are real 
values and, in general, are unequal. The main task of this 
paper is to investigate how the distributions of d™ and d""* 
will affect the network synchronization. 

We start by considering an unweighted, symmetrical net- 
work described by adjacency A = {anm,}, with Unm — 1 
if nodes n and m are connected, a„m = otherwise, and 



with 



0. The degree of node n is kn = Yl 
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m=l "n™ 



. To 

introduce the coupling gradient, we transform the adjacency 
matrix A as follows. For each pair of connected nodes n and 
m in the network, we deduce an amount g from a„m (the cou- 
pling that n receives from m) and add it to a„i„ (the cou- 
pling that m receives from n). In doing this, the total coupling 
strength between n and m is keeping unchanged. However, a 
coupling gradient will be generated, which is pointing from n 
to m. Denoting the resulted matrix as S. The coupling matrix 
C is then defined as c„m = knSnm/ J2i=i ^nj- The coupling 
gradient from n to m therefore is Acm„ = Cmn — c„m ~ 

f^mSm7i/ z-^j—1 ^'nij i^nSmn/ z-^j—i ^^j- ^^ practice, tiie di- 
rection of the coupling gradient is determined by a scalar field 
which, in terms of synchronization, is usually adopted as the 
node degree Jj, S '10]. Without losing the generality, we 
arrange coupling gradient flow from the larger-degree to the 
smaller-degree nodes on each link (the inverse situation can 
be achieved by setting g < 0). 

Now we discuss how the change of the gradient parameter g 
will affect the network synchronization. Noticing that in Eq. 
(O the value of ai is independent of g and, in constructing 
matrix C, we have (d™) = (fc) which is also independent 
of g. Therefore the introduction of gradient will only affect 
the value of ^(i™d°"*) in Eq. ([3]). Rearranging node index by 
a descending order of node degree, i.e. ki > k2 ■ ■ ■ > few, 
then the outgoing couplings of node n can be divided into 
two groups. Neighbors of node index m < n have elements 
Smn = 1 ~ 5 in matrix S and Cmn < 1 in matrix C (coupling 
gradient points to n), while for nodes of index m > nwe have 
Smn = 1 + 5 in matrix S and c,„„ > 1 in matrix C (coupling 
gradient departs from n). By this partition, the total outgoing 
coupling strength d°"* of node n reads 



JOllt _ J, 



l-.g 



p,.. 



1+9 

n,. 



p.. 



(4) 



with Pi<„ (Pi>n) is the probability for a randomly chosen 
node to have degree larger (smaller) than that of node n, and 
^i = -^J2i=i Sij is the normalizing factor defined on node. 
In calculating fli, again, we can divide the neighbors of i into 
two groups. Nodes of index j < i have s^ = 1+ g and nodes 
of index j > i have Si,j — 1 ~ g. Based on this partition, we 
obtain 

a = l + .g(Pj<. -Pj>0, (5) 

with P the same meaning as that of Eq. (|4|. For heteroge- 
neous networks of degree distribution P(/c) — C/c^'', wehave 

^i = l + :frr {'^^i~'' ~ ^nral - ^mTn K with kmax and kjnm 

denote the largest and smallest node degrees of the network, 
respectively. Inserting this into Eq. (|4|i, we obtain 

d7' =kn[F + G„] (6) 
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Finally we have 



[F + G{k)]k^P{k)dk 



(8) 



(9) 



Eq. (|9]l is our main result which tells how the network syn- 
chronization (Ec) changes with the coupling gradient (g) and 
the network topology (7). 

From Eq. Q we know that the introduction of coupling 
gradient changes only the weights H — F + G of the out- 
going couplings at each node, while the total coupling cost 
of the network is keeping unchanged. That is to say, gradi- 
ent changes the distribution of H from an even form (H — 1 
in unweighted network) to an uneven form (H = H{g, k) 
in weighted network). Physically, we can regard F as the 
symmetrical part of the couplings on each link, i.e. F ^ 
min(c„i, Cin), which depends only on the gradient parameter 
g and is decreased as g increases. In contrast, the term G is a 
joint function of g and A:„. While G increases with g monoton- 
ically, its exact value, however, are strongly modified by the 
node degree: node of larger degree assumes larger G [Eq. ([8]l]. 
According to the value of H, we are able to divide the network 
nodes into two groups. Nodes which have degree larger than 
kc have H > I, while those have degree smaller than kc have 
H < 1. The critical degree kc can be calculated by requiring 
H = 1. Under the assumption of k^ax >> kmin, we get 
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The uneven distribution of H becomes even clear when con- 
sidering the extreme cases of fc w kmax and k « kmm- 
From Eq. ^ we have i?fe«fe„„, = ^{l+g)ln^ 



2g 



(l-3)lni±f. Clearly, H, 
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Since the sum of H is constant for the network. 



and Hkak 

i.e. X]i=i Hi ^ ^' '^hs gradient effect thus can be understood 
as a shifting of partial of the weight H from smaller-degree to 
higher-degree nodes. 

For scale-free networks generated by the standard BA 
growth model fT], we have k^ax ~ kminN^i~^ . Inserting 
this into Eq. ^ we obtain 



H = F -\n 



l-g + '29 



1-7' 



(11) 



which basically tells the following: for a fixed gradient 
strength g, increase of the network homogeneity, i.e. the de- 
gree exponent 7, will make the distribution of H more ho- 
mogeneous and, as a result, suppress the network synchro- 
nization. In other words, the value of £c increases with the 
increase of 7. 
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FIG. 1: (Color online) For scale-free network of 1500 nodes, average 
degree 400, and degree exponent 7 = 3, the variation of the squared 
order parameter r as a function of the coupling strength e in the 
region of e G [O.SEc, I.SEc] under gradient strengths g = 0.5 (the 
left symbol curve), (7 = (the middle symbol curve), and g — —0.5 
(the right symbol curve). Apparently, the onset of synchronization is 
shifted to small coupling strengths at larger values of g. Each data 
is averaged over 10 network realizations. The three line curves are 
plotted according to Eq. (O, which predicts the behavior of r^ rea- 
sonably well in the region of e £ [sc, O.SEc]- Specially, the numeri- 
cal results of the critical coupling strength Sc are in good agreements 
with the theoretical results predicted by Eq. ([3}. 



The effects of coupling gradient and network topology on 
network synchronization can be summarized as follows. The 
change of the gradient strength g or the degree exponent 7 
does not change the total coupling cost of the network, it only 
redistributes the weight of the outgoing couplings at each node 
according to its degree information. When gradient g > 
is introduced, the outgoing couplings of small-degree nodes 
(having degree k < kc) will be reduced by an amount and 
added to those of large-degree nodes (having degree k > kc). 
This will induce a heterogeneous distribution in H which in 
turn will decrease the value of Ec (see Eq. (O). This enhance- 
ment of network synchronization, however, is modulated by 
the network topology. By increasing the degree exponent 7, 
the distribution of H tends to be homogeneous (i.e. H ^ 1) 
and, consequently, the network synchronization is suppressed. 
These are the mechanisms that govern the functions of cou- 
pling gradient and network topology. Therefore our findings 
are: 1) by coupling gradient, synchronization in non-identical 
networks can be enhanced; and 2) in comparison with ho- 
mogeneous networks, heterogeneous networks take more ad- 
vantages from coupling gradient and, under the same gradient 
strength, are more synchronizable. 

We now provide the numerical results. The networks are 
generated by a generalized BA model [12], which is able to 
generate networks of varying degree exponent 7. The fre- 
quency distribution is given by p{ijj) — (3/4)(l — w^) for 
— 1 < w < 1 and p{ijj) — otherwise. The initial phase 
6 of each oscillator is randomly chosen within range [0, 27r]. 



FIG. 2: For a scale-free network of 5000 nodes, average degree 100, 
and degree exponent 7 = 3, the variation of the critical coupling Ec 
as a function of the gradient parameter g. The solid line represents 
the theoretical results predicted by Eq. ^. 



A transition time T = 100 is discarded, and the value of r^ 
is calculated over another period of T = 100. To show the 
gradient effects on network synchronization, we have first cal- 
culated the variations of the squared order parameter r^ as a 
function of the coupUng strength e under three gradient pa- 
rameters: g = 0, 0.5 and —0.5. (According to our definition, 
(7 < means that coupling gradient is pointing from small- 
degree to larger-degree nodes.) The results are plotted in Fig. 
[T] Clearly, Ec becomes smaller at larger g. The three lines in 
Fig. [T]are plotted by Eq. (|2]i, which fit well with the numeri- 
cal results in the neighboring region of the onset point. More 
importantly, the value of Sc is predicted precisely by Eq. (O. 
(The precision of the predications dependent on the system 
size and the coupling density, larger and denser networks give 
better results.) 

To have a global picture on the gradient effect, we plot Fig. 
|2] the variation of the critical coupling strength Ec as a func- 
tion of the gradient parameter g. It is shown that, as g changes 
from — 1 to 1, the value of £c is monotonically decreased. This 
process of synchronization enhancement is well captured by 
Eq. ([8]) (the line curve in Fig. |2]), especially when g « 0. 
Since in our analysis we have assumed the network to be very 
large and dense, the mismatch between the theoretical and nu- 
merical results in Fig. |2]is understandable. 

Simulations have been also conducted on the dependence 
of Ec on 7. By the generalized BA model flZ|, we vary the 
degree exponent 7 continuously from 3 to 25, while keeping 
the size and average degree unchanged. As the prediction of 
Eq. ( fTTT ). it is found that, for each value of g, the critical cou- 
pling strength Ec increases monotonically with the increase 
of the degree exponent 7 [Fig. [3]]. Specially, for the case 
of g = 1 X 10~^, the numerical results are in good agree- 
ments with the theoretical approximations of Eq. ( fTOl ). As 
g increases, the mismatch between theoretical and numerical 
results will be enlarged which, again, can be alleviated by in- 
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FIG. 3: (Color online) For a scale-free network of 5000 nodes, av- 
erage degree 100, the variation of the critical coupling strength Ec 
as a function of the degree exponent 7 under gradient parameters 
g — 1 X 10^^ (the upper symbol curve) and g — 0.2 (the lower 
symbol curve). In both cases, Ec increases with 7. The solid lines are 
the theoretical results predicted by Eq. dllb . 



creasing the network size and density. 

A few remarks are in order. Firstly, while our theory is 
derived from large and densely connected networks, the gen- 
eral finding that synchronization can be enhanced by gradi- 
ent couplings applies to any network. The interesting thing 
is that, for a network of given degree distribution (not lim- 
ited to the scale-free type), our theory predicts quantitatively 



how large an improvement the network could benefit for a 
given gradient. From the findings, we are able to not only 
point out the optimal configuration for network synchroniza- 
tion, which happens at g = I [as shown in Fig. |2l, but also 
have a systematic understanding on the transition of the sys- 
tem performance as a function of gradient, and, more impor- 
tantly, the underlying mechanisms that govern this transition. 
Secondly, although similar findings about gradient effects had 
been discovered previously in studying global synchroniza- 
tion in identical networks IsLlgL llOtl . our analysises, however, 
are focusing on the onset of synchronization in non-identical 
networks. Another difference is that, by adopting the general- 
ized Kuramoto model, we are able to show analytically how 
the change of the coupling gradient will affect the synchro- 
nization (Eq. ^ and Fig. |2]) and what is the role of network 
topology in this process ((Eq. ( fTTT ) and Fig. O. It is noticed 
that in Ref. [11] the authors found numerically that the onset 
of synchronization in scale-free networks happens in advance 
to that of homogeneous networks, which can be understood 
readily from Eq. ( fTTT i. 



In summary, we have studied the effects of coupling gradi- 
ent on the onset of synchronization in nonidentical complex 
networks and found that, by coupling gradient, the network 
synchronization can be significantly enhanced and, in com- 
parison with homogeneous networks, heterogeneous networks 
are more synchronizable. We hope these findings could give 
insights to the collective behaviors in realistic systems. 
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